JOURNAL OF RESEARCH of the National Bureau of Standards— A. Physics and Chemistry 
Vol. 69 A, No. 4, July-August 1965 

Correlation of Successive Atomic Steps in Crystals by 

Relaxation Mode Analysis 

A. D. Franklin* 

(March 5, 1965) 



The relaxation mode analysis is used to compute the average value of the cosine of the 
angle between successive steps of an ion diffusing by a vacancy mechanism on a lattice. 
The technique is applicable to tracers (self-diffusion) or impurities, on any lattice. The 
result, given in terms of the eigenvalues and eigenvectors of a secular equation, is 
suitable for machine computation. Sample computations for self -diffusion on an fee lattice 
illustrate the method. It is necessary, in order to do the calculations, to have a lattice 
of finite size, the boundary of which is chosen either totally reflecting or completely trans- 
parent to vacancies. For self-diffusion, the boundary can be centered on the vacancy, and 
convergence of the computed value with boundary radius is rapid. For impurity diffu- 
sion, the boundary must be centered on the impurity, and convergence is much slower. 



1. Introduction 

Diffusion of atoms or ions through a crystal by the 
vacancy or interstitialcy mechanisms proceeds not 
by random- walk steps but by steps that are to some 
extent correlated with each other [lj. 1 The degree 
of correlat ion bet ween successive steps can be repre- 
sented by cos a , where 6 a is the angle between the 
paths of the ath step and the (o; + l)th step. If the 
term " tracer " is defined to mean any moving particle 
in the crystal, then where the tracer is a host ion 
movin g upon its own sublattice, as in self-diffusion, 
cos 6 a depends only upon the geomet ry of the lattice. 
If the tracer is an impurity, cos 6 a depends as well 
upon the way in which the impurity influences the 
motion of the vacancy in its vicinity. 

A numb er of m ethods have been developed for 
computing cos 6 a [2-6]. As shown by Compaan and 
Haven [3], the computations are very much simplified 
if use is made of the crystal symmetry. It is the 
purpose of this paper to show how maximum use 
can be made of the symmetry by employing the 
technique of rel axation mode analysis [7-10] to 
calculate cos 6 a . The method given is perfectly 
general, applying to both self-diffusion and impurity 
diffusion in any crystal lattice by the vacancy 
mechanism. With minor modification it could be 
applied to the interstitialcy mechanism as well. It 
is cast in a form that is extremely convenient for 
machine computation. 

*Work done while on a training assignment to the Theoretical Physics Division, 
UKAERE, Harwell, England. Permanent address: National Bureau of 
Standards, Washington, D.C. 

i Figures in brackets indicate the literature references at the end of this paper. 



2. Theory 

The development of the method will be given for 
any lattice, but the fee lattice will be used to illustrate 
the development and the method itself. The prob- 
lem can be visualized with the aid of figure 1. This 
shows the tracer at the center of a cube, with the 
nearest neighbor sites in the fee lattice in the middle 




Figure 1. Successive tracer steps in the fee lattice. 

The tracer is at the center, after the ath step which has left the vacancy at P. A 
possible (a+1) step is shown, to one of the twelve nearest neighbor sites , making 
an angle Oak with the ath step. For self-diffusion, only two kinds of nearest 
neighbor sites, here marked 1 and 2, need be considered. (Eq 12.) 
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of the edges of the cube. The tracer has just com- 
pleted the step shown by the arrow labeled a, leaving 
the vacancy in the position labeled p. In addition, 
one of the twelve possibilities for the (a+l)th step 
is shown. After the aih step, the vacancy will move 
in some fashion. It may immediately reexchange 
with the tracer, or it may first move directly or 
indirectly to one of the eleven other positions next 
to the tracer and then reexchange with the tracer. 
Or it may move away altogether and not return, so 
that the next step of the tracer is an exchange with 
a different vacancy. A different vacancy will 
approach the tracer in a random manner, provided 
the original vacancy is far away when the new one 
approaches, so that successive steps in which different 
vacancies figure are not correlated. Only steps 
involving t he sam e tracer need be considered in 
calculating cos d a for low vacancy concentrations. 
The problem is solved when the probability of each 
one of the twelve candidates for the (a + l)th step 
is calculated. Although figure 1 shows only one set 
of lattice sites from which the vacancy can make a 
direct exchange with the tracer, in general there may 
be several, the sites in each set being crystallographi- 
cally equivalent in the presence of the tracer. The 
solution of the problem can be found either by con- 
sidering the average behavior of a large number of 
identical systems, or by thinking in terms of the 
probability of finding the vacancy on a given site 
for one system, as done by Bardeen and Herring [1] 
and others. This latter view will be taken here. 

Let the probability that the vacancy occupies the 
kth site at a time t after the aih step be P k (t). Also, 
let w k0 be the frequency with which a vacancy oc- 
cupying the kth site will change places with the 
tracer at the origin. Then 



cos 6 a 



dak L 



=S w*o cos Bak I P k (t)dt. 



(la) 



In the relaxation mode analysis, eq (la) can be 
written as a vector equation, involving vectors in 
site-space, to be defined below. One vector contains 
all of the time-dependence of the probabilities P k (t), 
the other the geometrical quantities cos 6 ak . When 
the basis for these vectors is properly chosen, the 
time dependence becomes a simple series of ex- 
ponential decay functions, for which the integrations 
are readily performed, and eq (la) reduces to a 
rather simple form. 

Site space may be defined by considering each site 
to correspond to a vector in a basis set that spans 
the space. Any set of scalars assigned to the sites 
can be represented by a vector in this space, with 
the scalars for coordinates. If there are N accessible 
lattice sites, then the site space is JV-dimensional. 
Equation (la) may then be written rathei com- 
pactly in terms of two vectors 



cos 6 a = 



j; 



■ P(t)dt, 



(lb) 



where the coordinates of r' are w k0 cos 6 a jc and the 



coordinates of P(t) are the time-dependent probabili- 
ties of finding the vacancy on the various sites. 

The process of identifying the proper basis can be 
started by dividing the lattice sites into shells, the 
sites in each shell being equivalent to each other 
under the site-group, which is the group of rotational 
(proper and improper) symmetry operations of the 
crystal about the tracer as origin. Each shell defines 
a subspace of the whole site-space, each subspace 
being orthogonal to all the subspaces defined by the 
other shells. A vacancy on any site in one of the 
first few shells, counting outward from the origin, 
can make a direct exchange with the tracer. For a 
vacancy on a site in the nth such shell, the frequency 
w n0 with which the vacancy would change places with 
the tracer is the same for all sites. Hence r can be 



written 



= S W n< 



T n 



where the summation is over all shells from which 
direct exchange of the vacancy with the tracer can 
occur, and r' n is the vector whose coordinates are 
cos dak for all sites in the nth shell and zero for all 
others. Then eq (lb) becomes 



COS 0« = ]C 1V 7 



•f 



P(t)dt. 



(lc) 



The vector P(t). Since all of the time-dependence 
is contained in the vector P(t), only this need be con- 
sidered in identifying the desired basis. The com- 
ponents of P(t), which are the probabilities P k (t), are 
governed by the set of N rate equations 



^=-2 {w kl P k ~w lk P l ) 



(2) 



where w kt is the frequency with which a vacancy on 
the kth site will make a single, direct move to the Zth 
site. In this problem, the tracer at the origin is 
treated as a sink for probability — the probability 
of a vacancy moving to the origin from an adjacent 
site is positive, but the probability of the return step 
is taken as zero. Equation (2) can be written in 
terms of the vector P(t) in site space: 



dP 
dt 



= CP(t) 



(3) 



where (7 is the matrix of the w k i in eq (2). The solu- 
tions to eq (3) are 

P(0=S S /*u w exp (-\ gs t) (4) 

q s 

where the j qs are defined by the initial distribution 
of probability, at t = 0, the u qs are the eigenvectors 
of the matrix (7, and the \ 9 are the corresponding- 
eigenvalues, which are also the inverse relaxation 
times. The subscripts q and s are defined below. 

Use is made of the irreducible representations of 
the site-group in finding the eigenvectors and eigen- 
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values. The subscript s denotes the various irre- 
ducible representations, some of which may be of 
dimension greater than one. For these, several 
values of s, one for each dimension, will be required. 
The subscript q denotes one of the several eigenvec- 
tors (and eigenvalues) belonging to one dimension of 
a given irreducible representation. To each dimen- 
sion of an irreducible representation of the site- 
group, there corresponds a subspace of the whole 
site-space, spanned by a set of basis vectors \ m8 . 
Each eigenvector u, can be written as a linear 
combination of the basis vectors in one set: 



" / i drnqs * ms* 



(5) 



Any basis vector \ ms lies wholly within the subspace 
defined by one shell of sites, although several basis 
vectors in a set may be drawn from the subspace of 
the same shell. In this case, the shell has been 
divided into subshells by the irreducible representa- 
tion, and the indices n in eq (lc) and m in eq (5) 
specify the subshell. 

The constants/^, etc., can be found from the initial 
position x , y , Zo of the vacancy. At t = 

P(o)=SS/< s u ?s 

Q 8 

so that, since the eigenvectors are all mutually 
orthogonal, and using eq (5), 

f, s . = P(0) • u (JS =2 <W P(0) • V m5 . 



But P(0) is a vector in site-space with unity assigned 
to the initial site of the vacancy, and zeros elsewhere. 
If this initial site is labeled the kth site in the pth 
subshell, then P(0) is orthogonal to all of the basis 
vectors except the Y ps , and 



Jqs — ^pqs\ V ps)k 



(6) 



where (V ps ) k is the kth component of the basis vector 
\ ps . Combining (4), (5), and (6) the probability 
vector P(t) may be written 

P(0=S 2D X) a pqg a mq8 (V P8 ) k V ms exp (— \ qst ). (7) 

q s m 

To recapitulate the subscripts, s specifies the irreduc- 
ible representation, or specifies one of the dimensions 
of a multidimensional irreducible representation, and 
m both denotes a basis vector within the set belong- 
ing to the sth irreducible representation and also 
identifies a subshell of lattice sites. Together, these 
subscripts cover the whole dimensionality of the site 
space. The subscript q indicates one of the eigen- 
vectors belonging to the sth representation — it runs 
parallel with but is independent of m. The initial 
position of the vacancy is indicated by p (one of the 
set of subshells identified by the running index m) 
and k, which gives the site within the subshell. 

The vector r n . The vectors r' n may also be 
expressed in terms of the same basis. Each set of 



basis vectors reflects the symmetry associated with 
the corresponding irreducible representation. In 
particular, there will always exist 3 sets, reflecting x, 
//, and z symmetry in real space, although they may 
not necessarily all correspond to the same irreducible 
representation as they do for cubic crystals. Any 
function with vector symmetry in real space will, 
when evaluated at each lattice site, generate a vector 
in site space that can be expressed as a linear com- 
bination of the basis vectors drawn only from the x, 
y, and z sets. For these vectors the subscript s will 
be replaced by x, y, or z. 

Using column vectors in site space, the x, y, and z 
basis vectors for the nth subshell may be written 



3m 



etc. 



where the x nk , y nk , and z nk are the coordinates in real 
space of the kth site in the nth subshell, and the L nx , 

etc, are the normalization constants ( T^, xl* J > 

etc. Now the vectors r n in site space are generated 
by the function cos 9 ak evaluated at the lattice sites 
of the nth subshell in real space. This function can 
also be written 



r<)-r, 



COS ak = 

\ l 'o\ |r*| 

where r is the vector in real space with coordinates 
x , y 0y and z from the origin to the site from which 
the tracer made the ath step (the initial site of the 
vacancy at t=0), and r^ is the vector to the kth 
site. Then the vector r n} in site space, can be 
written 



1 


1 


Tn—~ 


kolknl 



{ XoL nx \ nx -\-yoL n y\ n y-\-z G L NZ \ f/z 



(8) 



Cos $ a . With eq (7) for P(t) and eq (8) for r' n , 
the scalar product of these two vectors can be written 

r' n -P(t)= -lL{x 2(I « A " 6XP ( - W) 



\r \ 



J px g 



j x i Qjpqy&nqy ^P ^ h ny t) 



Zo-Lnz 



2 (Lpq/lnqz eX P ( 



~\ut) y 



where use has been made of the mutual orthogonality 
of the basis vectors V m . s . When this is inserted into 
eq (lc) for cos 6 a and the integration performed, 
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the result is 



COS0 a = 




0) 



This equation is in a form that is quite suitable for 
machine computation. To the irreducible repre- 
sentations denoted by x, y, and z, there correspond 
secular equations which can be derived by standard 
methods. The qih solution to these secular equa- 
tions provides values for \ gx , etc. Each solution 
also gives a set of coefficients, a nqx , etc., which are 
the contributions from the basis vectors \ nx to the 
gth eigenvector u gx . The vacancy is initially in the 
pih subshell with coordinates x Q , y , and z in real 
space. The summation denoted by n is taken over 
all subshells from which the vacancy can make 
direct exchanges with the tracer, and by q over all 
solutions to each secular equation. 

Equation (9) may be used to obtain the quantities 
appearing in Mullen's [5] correlation functions for 
diffusion in anisotropic lattices. Mullen showed 
that the important quantities are of the form 



*xa,a)- 



A A 

Jy ««// a-\-l 



«?(A B) = 



(xi) 2 



where xi is the projection of the ath tracer step on 
the x-principal diffusion axis, and the superscript A 
indicates that this step is of the A type. The dis- 
tinction between the A, B, etc., types of steps is the 
same as the distinction among the subshells denoted 
in eq (9) by the index n. From eq (9) it follows that 



e?(A, B)- 



Wbo (x A ) 2 L Bx 



r A r B 



Q 



&Aqx&Bqx 



where A, B as subscripts denote the subshells corre- 
sponding to the A, B types of step. 

Boundary conditions. In performing the computa- 
tions, it is necessary to restrict the vacancy to a 
finite number of sites around the tracer. This im- 
poses a boundary and necessitates a choice of bound- 
ary conditions. A convenient choice is to take a 
symmetrical boundary by including all sites out to a 
certain number of steps from the origin. Possible 
jumps of the vacancy to sites beyond the boundary 
can then be treated in one of two simple ways. If 
these jumps are allowed, but the return jumps are 
not, then the boundary is " transparent" and acts as 
a perfect sink. This is the boundary condition used 
by Compaan and Haven [3] and by Mullen [5]. Use 
of the transparent boundary neglects reflection of the 
vacancy from the region beyond the boundary. 
Since the path of the reflected vacancy would, to 
some extent, be correlated with the initial step, the 



transparent boundary condition underestimates the 
value of cos 6 a . On the other hand, the jumps 
through the boundary can be omitted altogether from 
the matrix (7, in which case the boundar y is tot ally 
"reflecting" and the result overestimates cos 6 a . If 
the computations are performed for both conditions, 
then an upper and a lower limit are obtained for 
cos 6 a . 

3. Self-Diffusion in the FCC Lattice 



The method can be illustrated by calculating cos 6 a 
for self- diffusion on a fee lattice. Compaan and 
Haven [3], using a somewhat different technique to 
solve essentially the same equations, have given a 
value for cos 6 a of —0.12268 for this case, to which 
the results of the present method may be compared. 

For the cubic lattices, the x, y, and z irreducible 
representations form a degenerate set, with the same 
secular equations. For illustrative purposes, we may 
consider the particular case of the fee lattice, and 
restrict direct exchanges to the nearest neighbor shell 
only, without subshells. Then p=n=l, r =r n 
= (xl+yl-{-ziy /2 7 and eq (9) reduces to the very 
simple form 

(10) 



cos# a =-w 10 y^ 



Furthermore, a simple prescription can be given for 
the secular equation for the x, y, z irreducible repre- 
sentation (T lu in the notation of Eyring, Walter, and 
Kimball [12]) in the fee lattice. This equation is 

C'u q =\ Q u Q 

where C is a matrix whose entries are: 

c^ = sum of the jump frequencies of all possible jumps 
of the vacancy out of a site included in the ith 
basis vector, into sites not included in the same 
vector. 
Ci j = —gij(d j /d i y /2 w ij , 

where g tj is the number of sites included in the 
ith basis vector from which jumps into any 
site included in the jth basis vector are made; 
dj and d t are the number of nonzero entries in 
the ^'th and ith basis vectors, respectively; and 
Wtj is the jump frequency for jumps of the 
vacancy from a site in the iih basis vector 
to one in the ^'th basis vector. 

In performing the calculation, the boundary is 
chosen, thus limiting the number of lattice sites 
considered explicitly and therefore the dimension of 
the matrix C'. The prescription given above is 
used to write down the matrix, and its eigenvalues 
and eigenvectors found by numerical computation. 
The eigenvalues \ q and the coefficients a iQ of the first 
basis vector (corresponding to the sites that are 
nearest neighbors to the tracer) in the eigenvectors 
are inserted into eq (10) to give cos 6 a . 

There is a convergence problem in these calcula- 
tions, arising from the presence of a boundary of 
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finite radius. For self-diffusion, there is no interac- 
tion between the tracer and the vacancy to bind 
them together, so that the vacancy may be expected 
to make long excursions away from the tracer, and 
the presence of the boundary may influence the 
computed value of cos 6 a even for a boundary quite 
far from the origin. If the boundary be made 
spherical, an estimate of the rate of convergence of 
the calculated value of cos 6 a as the radius R of the 
boundary increases may be made. For this estimate, 
use is made of the analogy pointed out by Ham [13] 
between the solutions to a steady-state diffusion 
problem and the solutions to an electrostatic problem 
of similar geometry. In the form given by Com- 
paan and Haven [13] eq (2) can be integrated with 
respect to time from to oo and written 



r 



^f=i=wS(F,-F P ); 



/. 



o dt 



dt=0=w^(V l -V j ), j^p (11) 



where w is the single jump frequency needed in the 
self-diffusion case, the subscript p denotes the initial 
position of the vacancy, and 



v t = r Ptifidt. 



Compaan and Haven noted that eqs (11) are identical 
to the equations governing current flow in a resistor 
network. However, these equations apply equally 
to the charge and potential distribution in a capacitor 
network. Let the lattice sites be connected by 
capacitors of capacitance w, and let the potential at 
the ith site be V t . Then if a unit charge is placed 
on the pih site, and zero charge on all others, the 
potentials on the lattice sites will conform to eqs 
(11). For the fee lattice, as shown by Compaan and 
Haven, cos 6 a is given by 



cos 6 a =—w(V 1 +2V 2 ) 



(12) 



where Vi and V 2 are the potentials at the sites 
marked 1 and 2 in figure 1, when a positive unit 
charge is at site p and a negative one at site 1. 

An estimate of these potentials may be made when 
the radius of the boundary is large enough to allow 
replacing the discrete network by a homogeneous 
dielectric sphere, with a radius R and with a charge 
— q placed at ^2a and another one q placed at 
— V2a, where a is the cation-anion separation (fig. 2). 
The origin is fixed at zero potential, and the boundary 
also taken as fixed at zero potential to reproduce the 
transparent boundary condition described above. 
The problem is then that of the potential due to 
fixed charges within a conducting spherical shell held 
at zero potential. As shown by Abraham and 
Becker [14], the same potential would be produced 
by image charges ±Q=±qR/^J2a } located at dis- 




Figure 2. Electrostatic model of the transparent boundary. 

R is the radius of the boundary; r and </> define the position of an arbil rary point 
within the boundary. Real charges— q and q are placed at the initial vacancy 
site and the inversion-related position, a distance -\J2a from the origin. Image 
charges Q and — Q are placed on the 7 axis ±R 2 /^2a from the origin. 



tances ±R 2 /^!2a from the origin as shown in figure 
2. If R^>a, then the potential at an arbitrary 
point within the sphere specified by spherical coordi- 
nates r and is 



"="•+!{ 



8a 6 



cos <^>+^ (5 cos 3 0—3 cos </>)... V 

(13) 



in which the first term V is due to the real charges 
within the sphere, and the second term is due to the 
image charges, and therefore represents the contribu- 
tion of the boundary. 

When potentials at sites 1 and 2 are evaluated 
using eq (13) and inserted in eq (12), the result is 



cos 



0«=cos0 aO +constantX (jkj ~^{%) • • • ' ( 14 ) 



This equation may be used as an extrapolation 
fo r nmla t o obtain a value for cos d a0 from the values 
of cos 6 a computed for finite boundaries with eq (10). 
It suggests that the computed value, for large enough 
values of the radius R of the boundary sphere, 
should approach the true value as l//i >3 . Values of 
cos 6 a , computed for both reflecting and transparent 
boundary, are plotted in figure 3 against (R/a)~ 3 , 
where R has been identified with the distance from 
the origin to the closest sites from which return 
jumps were not allowed (transparent boundary) or 
to which forward jumps were not allowed (reflect ing 
boundary). The values for the transparent bound- 
ary show the expected i?~ 3 behavior for large R. If 
eq (14) is used as an extrapolation formula, with 
arbitrary co efficien ts, the value of —0.12266 is 
obtained for cos 6 aQ , in essential agreement with the 
value computed by Compaan and Haven [3], 

This R~ z convergence is independent of the specific 
lattice considered, and arises simply from the fact 
that the tracer position after the ath step is taken as 
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/ REFLECTING BOUNDARY 



* BOUNDARY CENTERED ON VACANCY 



TRANSPARENT BOUNDARY 



5 10 

(R/a)" 3 XIO 3 

R= RADIUS OF BOUNDARY 



Figure 3. Dependence of cos B a upon (R/a) -3 , for both trans- 
parent and reflecting boundaries centered on the tracer, and for 
a transparent boundary centered on the vacancy. 

The solid line for the transparent boundary was drawn using eq (14) as an extrap- 
olation formula, plotting only the term in R~K 



the origin, around which a symmetrical boundary 
is placed. For impurity diffusion in which the 
impurity ion and the vacancy are bound together, 
the excursions of the vacancy out to the boundary 
will be less likely than in the case of self-diffusion, 
with no binding. Hence the convergence for im- 
puritv diffusion should be more rapid than given 
by R~\ 

It is possible to achieve more rapid convergence 
in the special case of self-diffusion by placing the 
position of the vacancy after the cutri step at the 
origin, as done by Compaan and Haven [3] and by 
Mullen [5], although this cannot be done for impurity 
diffusion. For self-diffusion the tracer is exactly 
like all other atoms or ions in the crystal. Compaan 
and Haven showed that for crystals in which the 
tracer site is a center of symmetry or has a mirror 
plane, a two-fold, or a three-fold axis of rotation 
perpendicular to the direction of the ath step, the 
trajectories of the vacancy that cross the tracer need 
not be terminated at the tracer. In the electrostatic 
model, this is equivalent to removing the condition 
that the tracer site be held at zero potential, and in 
the language of the relaxation mode analysis, it is 
equivalent to including jumps of the vacancy out of, 
as well as into, the tracer site, so that the tracer no 
longer is a sink for the probability. Under these 
circumstances, the problem has the symmetry of 
the site-group at the initial position of the vacancy, 
and this position can be taken as the origin, with 



the boundary distributed symmetrically about it, 
rather than about the tracer. In the electrostatic 
model for the fee lattice, the spherical boundary 
can be taken about the initial position of the vacancy 
as a center, and only at this point is there any charge 
within the sphere. Such a boundary simply adds 
a con stant to the potential everywhere so that 
cos B a , which depends only upon differences in poten- 
tial, is not affected by the presence of the boundary 
in the continuum approximation . For a real lattice, 
there will be a dependence of cos 6 a for boundary 
radii too small to allow a continuum approximation. 
The rapid convergence that results when the 
boundary can be centered on the vacancy can be 
illustrated by direct calculation. In the relaxation 
mode analysis of this case, the initial distribution 
of probability is spherically symmetric, and the only 
modes excited will be the spherically symmetric ones, 
the A lg for the fee lattice. Equation (10) will be 
replaced by 



COS d a = W 23 






4 






+ 



2a< 



n, 



2a 4g 

V^4 



fl5g \ 



where a iq are the coefficients of the ith basis vector 
in the qili eigenvector of the A Xg irreducible repre- 
sentation, \ q are the corresponding eigenvalues, and 
Ui is the number of sites in the ith. shell (there are 
no subshells for the spherically symmetric modes). 
Using this equation, and setting a transparent bound- 
ary at R/a= -\/4:2 a nd R/a=Tf52, we obtained 
-0.12268 for cos d a in both cases (fig. 3). 



The 



value had converged even for the smaller boundary, 
and agrees exactly with Compaan and Haven's 
value. 

It should be noted that the boundary used in this 
calculation was not exactly identical to that used 
either by Compaan and Haven [3] or by Mullen [5]. 
In the present calculations, the boundary was taken 
as spherical. The distance from the origin to any 
point on the boundary was the same in all directions, 
but the number of steps required for the vacancy to 
reach the boundary was not the same. If, as done 
by Compaan and Haven, and Mullen, the boundary 
was taken as a certain number of steps from the 
origin in all directions, the secular equation could 
have been somewhat smaller in size and the problem 
somewhat simpler, without affecting the convergence 
significantly. 

The present method of calculating cos 6 a , using 
the relaxation mode analysis, is basically similar to 
Mullen's [5] method. Equation (4) of the present 
paper, giving the probability vector T?(t), can be 
made equivalent to Mullen's eq (10) by an integra- 
tion over t from zero to infinity and multiplication 
of each component by the sum of the frequencies 
of all jumps out of the corresponding site. Mullen 
treats only the case of self-diffusion, with the bound- 
ary centered on the vacancy. Simplification of the 
equations is brought about by introducing only sites 
not symmetry-related to each other under the 
crystallographic point group. This is equivalent to 
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the choice of the A ig mode as described above. The 
two treatments lead to matrices of the same dimen- 
sionality. In Mullen's formulation, it is necessary 
to invert the matrix; here, the eigenvalues and 
eigenvectors must be calculated. 

For impurity diffusion, the relaxation mode 
analysis simplifies the equations by expressing the 
probability vector P(t) as a sum of eigenvectors of 
the matrix C, each multiplied by its relaxation func- 
tion. Equation (lb) then selects only those eigen- 
vectors not orthogonal to r f . The same procedure 
could be applied to Mullen's eq (10). The initial 
probability distribution, r(o), would be expressed as 
a sum, each term of which would be equivalent to an 
eigenvector of Mullen's T matrix. The choice of 
the eigenvectors as a basis would reduce T to its 
irreducible representations, of which only the vector 
representation would have to be kept, thus achieving 
the same reduction in the size of the matrix to be 
handled as obtained in the relaxation mode analysis. 
In addition, the boundary would have to be centered 
on the site of the impurity, since this is the only 
center of symmetry in the presence of the impurity. 

4. Summary 

Correlation between successive steps of an impu- 
rity or tracer as it diffuses through a crystal by the 
vacancy mechanism can often be expressed in terms 
of cos 6 a , the average of the value of the cosine of 
the angle between the ath and (a-fl)th steps, 
cos 6 a can be computed for any lattice using the 
methods of relaxation mode analysis, resulting in 
eq (9) of this paper, and reducing the problem to a 
form well suited to machine computation. A sym- 
metrical boundary is taken about the impurity or 
tracer. When the boundary is transparent, so that 
the vacancy does not return once it has penetrated 
the boundary, th e comp utation gives a lower limit to 
the magnitude of cos 6 a . An upper limit is produced 
when the bounda ry is t otally reflecting. The com- 
puted values of cos d a for self-diffusion with the 
transparent boundary converge as 1//Z 3 , where R 
is the distance from the origin to sites just beyond 
the boundary. For impurity diffusion in which the 
vacancy and impurity are attracted to each other, 
the convergence should be somewhat more rapid. 



The method has been illustrated by computing 
cos 6 a for self -diffusion in the fee lattice. The 1/A >3 
convergence was demonstrated, and an extrapolated 
value of 0.12266 obtained, in essential agreement 
with the value of 0.12268 given by Compaan and 
Haven [3]. For the special case of self-diffusion, 
Compaan and Haven showed that for sufficiently 
symmetric crystals, the origin can be placed at the 
initial position of the vacancy and a symmetrical 
boundary placed around it, producing much more 
rapid convergence. Computations for this model 
were also made for the fee lattice, with two sizes of 
boundary, and both gave values in exact agreement 
with Compaan and Haven's value. 



The author is indebted to A. B. Lidiard for the 
original suggestion that the relaxation mode analysis 
could be used to calculate cos d ay and to R. C. F. 
McLatchie for performing the computations. 
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